xtset id year
xtreg Carbon did  i.year,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
xi:xtscc Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe vce(bootstrap,reps(1000))
****PSMDID*********
set seed 13579
gen tmp = runiform()
sort tmp
psmatch2 treat Ingy InFDI Inbus InGDP Inurba, out(Carbon) logit ate neighbor(4) common caliper(.05) ties
pstest Ingy InFDI Inbus InGDP Inurba, both graph  //最近邻1：4/失败
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if _weight !=.,fe r

psmatch2 treat Ingy InFDI Inbus InGDP Inurba, out(Carbon) logit radius ate  common caliper(.005) ties
pstest Ingy InFDI Inbus InGDP Inurba, both graph  //半径匹配
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if _weight !=.,fe r

psmatch2 treat Ingy InFDI Inbus InGDP Inurba, out(Carbon) logit kernel ate  common ties
pstest Ingy InFDI Inbus InGDP Inurba, both graph  //核匹配
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if _weight !=.,fe r
*****工具变量******
gen Inhlw=log( 互联网宽带接入用户数_万户_全市 )
gen hlw=Inhlw*time
xi:xtivreg2 Carbon Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba (did= hlw) i.year,fe first r
gen Inphone=log( 移动电话年末用户数_万户_全市 )
gen phone=Inphone*time
xi:xtivreg2 Carbon Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba (did= phone) i.year,fe first r  //民众高效的信息通讯要求政府同样提高效率
****平行趋势检验******
gen distance = year - first_treat
forvalues i=1/11 { 
gen pre_`i'  = 0 
replace pre_`i'  = 1 if treat== 1 & distance== -`i'
}
forvalues i=1/2 { 
gen after`i'  = 0 
replace after`i'  = 1 if treat== 1 & distance== `i'
}
gen current  = 0
replace current = 1 if treat== 1 & distance== 0
xtreg Carbon pre_5 pre_4 pre_3 pre_2 pre_1 current after1 after2 Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year, fe r
coefplot,keep(pre_5 pre_4 pre_3 pre_2 pre_1 current after1 after2 after3 ) levels(95) vertical lcolor(black) mcolor(black) msymbol(circle_hollow) ytitle(Coefficient, size(small))  ylabel(, labsize(small) angle(horizontal) nogrid) yline(0, lwidth(vthin)lpattern(solid) lcolor(black)) xtitle(Years since the policy implementation, size(small)) xlabel(,labsize(small))  graphregion(fcolor(white) lcolor(white) ifcolor(white) ilcolor(white)) ciopts(recast(rcap)) xline(10.5, lwidth(vthin) lpattern(solid)lcolor(black))
****didmultiplegt稳健性***
set seed 13579
did_multiplegt Carbon id year did , robust_dynamic dynamic(2) placebo(5) breps(50) cluster(id) controls( Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba)
event_plot e(estimates)#e(variances), default_look graph_opt(xtitle("Years since the policy implementation") ytitle("Average causal effect") ///
  title("Parallel trend test under DIDM") xlabel(-5(1)2)) stub_lag(Effect_#) stub_lead(Placebo_#) together
****随机处理组*************
***500次***
xtset id year
global xlist " Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba  "
xtreg Carbon did $xlist i.year,fe r
*真实系数 = -0.245
*真实t值  = -3.08
*真实p值 = 0.002
mat b = J(500,1,0) //* 系数矩阵
mat se = J(500,1,0) //* 标准误矩阵 
mat p = J(500,1,0) //* P值矩阵
forvalues i=1/500{
	use "C:\Users\86189\Desktop\互联网+政务服务\随机处理组\data2.dta", clear
	xtset id year  //面板数据声明
	keep if year==2013   //保留一期数据
	sample 80, count   //随机抽取80个城市
	keep id  //得到所抽取样本的id编号
	save match_id.dta, replace   //另存id编号数据
	merge m:m id using "data1.dta" //与原数据匹配
	gen new_treat= (_merge == 3) //将所抽取样本赋值为1，其余为0，得到政策分组虚拟变量
	gen new_time= (year >= 2017) //生成政策时间虚拟变量
	gen new_did = new_treat*new_time
	xtreg Carbon new_did $xlist i.year,fe r
	* 将回归结果赋值到对应矩阵的对应位置
	mat b[`i',1] = _b[new_did] //* 系数矩阵
	mat se[`i',1] = _se[new_did] //* 标准误矩阵
	* 计算P值并赋值于矩阵
	mat p[`i',1] = 2*ttail(e(df_r), abs(_b[new_did]/_se[new_did]))
}
* 矩阵转化为向量
svmat b, names(coef)
svmat se, names(se)
svmat p, names(pvalue)
drop if pvalue1 == .
label var pvalue1 p值
label var coef1 估计系数
keep coef1 se1 pvalue1  
save placebo.dta,replace

use placebo.dta,clear
twoway (kdensity coef1,yaxis(1)) (scatter pvalue1 coef1, yaxis(2) msymbol(smcircle_hollow) mcolor(blue)), ///
title("Placebo Test") ///
xlabel(-0.3(0.1)0.3) ylabel(,angle(0)) ///
xline(-0.245, lwidth(vthin) lp(shortdash)) xtitle("Coefficients") ///
yline(0.002,lwidth(vthin) lp(dash)) ytitle(p value) ///
legend(label(1 "kdensity of estimates") label( 2 "p value")) ///
plotregion(style(none)) ///无边框
graphregion(color(white)) //白底
***1000次***
xtset id year
global xlist " Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba  "
xtreg Carbon did $xlist i.year,fe r
*真实系数 = -0.245
*真实t值  = -3.08
*真实p值 = 0.002
mat b = J(1000,1,0) //* 系数矩阵
mat se = J(1000,1,0) //* 标准误矩阵 
mat p = J(1000,1,0) //* P值矩阵
forvalues i=1/1000{
	use "C:\Users\86189\Desktop\互联网+政务服务\随机处理组\data2.dta", clear
	xtset id year  //面板数据声明
	keep if year==2013   //保留一期数据
	sample 80, count   //随机抽取80个城市
	keep id  //得到所抽取样本的id编号
	save match_id.dta, replace   //另存id编号数据
	merge m:m id using "data1.dta" //与原数据匹配
	gen new_treat= (_merge == 3) //将所抽取样本赋值为1，其余为0，得到政策分组虚拟变量
	gen new_time= (year >= 2017) //生成政策时间虚拟变量
	gen new_did = new_treat*new_time
	xtreg Carbon new_did $xlist i.year,fe r
	* 将回归结果赋值到对应矩阵的对应位置
	mat b[`i',1] = _b[new_did] //* 系数矩阵
	mat se[`i',1] = _se[new_did] //* 标准误矩阵
	* 计算P值并赋值于矩阵
	mat p[`i',1] = 2*ttail(e(df_r), abs(_b[new_did]/_se[new_did]))
}
* 矩阵转化为向量
svmat b, names(coef)
svmat se, names(se)
svmat p, names(pvalue)
drop if pvalue1 == .
label var pvalue1 p值
label var coef1 估计系数
keep coef1 se1 pvalue1  
save placebo.dta,replace

use placebo.dta,clear
twoway (kdensity coef1,yaxis(1)) (scatter pvalue1 coef1, yaxis(2) msymbol(smcircle_hollow) mcolor(blue)), ///
title("Placebo Test") ///
xlabel(-0.3(0.1)0.3) ylabel(,angle(0)) ///
xline(-0.245, lwidth(vthin) lp(shortdash)) xtitle("Coefficients") ///
yline(0.002,lwidth(vthin) lp(dash)) ytitle(p value) ///
legend(label(1 "kdensity of estimates") label( 2 "p value")) ///
plotregion(style(none)) ///无边框
graphregion(color(white)) //白底
****随机政策时间**********************
***500次***
mat b = J(500,1,0)
mat se = J(500,1,0)
mat p = J(500,1,0)
forvalues i=1/500{
 use "C:\Users\86189\Desktop\互联网+政务服务\随机政策时间\data1.dta", clear
 xtset id year
 sample 298, count
 keep year
 mkmat year, matrix(sampleyear) 
 use "C:\Users\86189\Desktop\互联网+政务服务\随机政策时间\data1.dta",clear
 xtset id year
 gen new_did = 0
 foreach j of numlist 1/298{
  replace new_did = 1 if (id == `j' & year >= sampleyear[`j',1] )
 }
 qui xtreg Carbon new_did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
 mat b[`i',1] = _b[new_did]
 mat se[`i',1] = _se[new_did]
 mat p[`i',1] = 2*ttail(e(df_r),abs(_b[new_did]/_se[new_did]))
}
svmat b, names(coef)
svmat se, names(se)
svmat p, names(pvalue)
drop if pvalue1 == .
label var pvalue1 p值
label var coef1 估计系数
keep coef1 se1 pvalue1  
save placebo.dta,replace 

use placebo.dta,clear
twoway (kdensity coef1,yaxis(1)) (scatter pvalue1 coef1, yaxis(2) msymbol(smcircle_hollow) mcolor(blue)), ///
title("Placebo Test") ///
xlabel(-0.3(0.1)0.1) ylabel(,angle(0)) ///
xline(-0.245, lwidth(vthin) lp(shortdash)) xtitle("Coefficients") ///
yline(0.002,lwidth(vthin) lp(dash)) ytitle(p value) ///
legend(label(1 "kdensity of estimates") label( 2 "p value")) ///
plotregion(style(none)) ///无边框
graphregion(color(white)) //白底
***1000次***
mat b = J(1000,1,0)
mat se = J(1000,1,0)
mat p = J(1000,1,0)
forvalues i=1/1000{
 use "C:\Users\86189\Desktop\互联网+政务服务\随机政策时间\data1.dta", clear
 xtset id year
 sample 298, count
 keep year
 mkmat year, matrix(sampleyear) 
 use "C:\Users\86189\Desktop\互联网+政务服务\随机政策时间\data1.dta",clear
 xtset id year
 gen new_did = 0
 foreach j of numlist 1/298{
  replace new_did = 1 if (id == `j' & year >= sampleyear[`j',1] )
 }
 qui xtreg Carbon new_did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
 mat b[`i',1] = _b[new_did]
 mat se[`i',1] = _se[new_did]
 mat p[`i',1] = 2*ttail(e(df_r),abs(_b[new_did]/_se[new_did]))
}
svmat b, names(coef)
svmat se, names(se)
svmat p, names(pvalue)
drop if pvalue1 == .
label var pvalue1 p值
label var coef1 估计系数
keep coef1 se1 pvalue1  
save placebo.dta,replace 

use placebo.dta,clear
twoway (kdensity coef1,yaxis(1)) (scatter pvalue1 coef1, yaxis(2) msymbol(smcircle_hollow) mcolor(blue)), ///
title("Placebo Test") ///
xlabel(-0.3(0.1)0.1) ylabel(,angle(0)) ///
xline(-0.245, lwidth(vthin) lp(shortdash)) xtitle("Coefficients") ///
yline(0.002,lwidth(vthin) lp(dash)) ytitle(p value) ///
legend(label(1 "kdensity of estimates") label( 2 "p value")) ///
plotregion(style(none)) ///无边框
graphregion(color(white)) //白底
*****固定效应及聚类*******
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id)
estimates store fe1
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year provinceid ) vce(cluster id)
estimates store fe2
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year provinceid#year ) vce(cluster id)
estimates store fe3
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year provinceid provinceid#year ) vce(cluster id)
estimates store fe4

reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id#year)
estimates store fe5
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id id#year)
estimates store fe6
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id provinceid year)
estimates store fe7
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id provinceid id#year)
estimates store fe8
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id provinceid year#provinceid)
estimates store fe9
reghdfe Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba ,a(id year) vce(cluster id provinceid year#provinceid id#year)
estimates store fe10
coefplot fe1 fe2 fe3 fe4 fe5 fe6 fe7 fe8 fe9 fe10 ,keep(did) xline(0) levels(99)  msymbol(D)
*****排除其他政策干扰*****低碳城市/碳市场
xtreg Carbon did LCC Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
xtreg Carbon did ETS Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
*******替换被解释变量*****
gen CI=log( carbon百万吨 / 地区生产总值_当年价格_亿元_市辖区 )
xtreg CI did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
gen PCE=log( carbon百万吨 / 年末户籍人口_万人_全市 )
xtreg PCE did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
*******异质性********
****东中西部*****
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if geo==1,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if geo==2,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if geo==3,fe r
*****城市规模*****
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if p==5,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if p==3,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if p==1,fe r
****经济发展程度******
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if g==5,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if g==3,fe r
xtreg Carbon did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year if g==1,fe r
****传导机制*****
****数字经济发展****
xtreg 数字经济水平_主成分 did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
xtreg 信息传输、计算机服务和软件业从业人员占比 did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
//或者//从业人口可作为信息产业基础
gen employee_ICT=log(信息传输、计算机服务和软件业从业人员占比)
xtreg employee_ICT did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
xtreg Structure did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
//或者//
xtreg 数字普惠金融指数 did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
xtreg Dig_eco did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r
***企业数字化转型**企业年报词频依地区取均值
xtreg dig did Ingy InFDI InFIN Insci InFin Inbus InGDP Inurba i.year,fe r


